
###################################################
#
#               Code to replicate 
#  The Influence Gap: Unequal policy responsiveness 
#              to men and women 
#   
#    Ruben Berge Mathisen (University of Bergen)
#
#       Tables and Figures for Appendix
#                April 2023
#
#   Original analyses conducted in R version 
#             4.1.3 (2022-03-10)
#
##################################################



#####################################################
######## Packages and data prep #####################
#####################################################

library(readxl)
library(dplyr)
library(ggplot2)
library(ggeffects)
library(interplot)
library(gridExtra)
library(plyr)
library(kableExtra)
library(stargazer)
library(interflex)
library(tidyr)
library(car)
library(RColorBrewer)
library(scales)
library(countrycode)
library(haven)

logitTransform <- function(p) { log(p/(1-p)) }

# Main datasets 
d_us <- read_excel("Main dataset for United States, from Gilens 2012.xlsx")
d_us$country <- "United States"

#dr <- read_excel("OpinionPolicyDataset FINAL v5 GENDER.xlsx")
d_no <- read_excel("Main dataset for Norway.xlsx")
d_no$country <- "Norway"

# Data on women in parliament
parl <- read_excel("Share of women MPs norway and us.xlsx")
## Fill in data for non-election years, Norway
parl_nor <- rbind.fill(parl%>%filter(country=="norway"),
                       parl%>%filter(country=="norway")%>%mutate(year=year+1),
                       parl%>%filter(country=="norway")%>%mutate(year=year+2),
                       parl%>%filter(country=="norway")%>%mutate(year=year+3))
## Make lagged versions of variable and 4-year averages, Norway
parl_nor <- parl_nor %>% 
  arrange(year) %>%
  mutate(share_t1=lead(share,1),
         share_t2=lead(share,2),
         share_t3=lead(share,3),
         share_t4=lead(share,4),
         share_avg4=(share+share_t1+share_t2+share_t3+share_t4)/5)
## Fill in data for non-election years, US
parl_us <- rbind.fill(parl%>%filter(country=="us"),
                      parl%>%filter(country=="us")%>%mutate(year=year+1))
## Make lagged versions of variable and 4-year averages, US
parl_us <- parl_us %>% 
  arrange(year) %>%
  mutate(share_t1=lead(share,1),
         share_t2=lead(share,2),
         share_t3=lead(share,3),
         share_t4=lead(share,4),
         share_avg4=(share+share_t1+share_t2+share_t3+share_t4)/5)
parl_full <- rbind(parl_nor,parl_us)
parl_full <- parl_full %>% 
  dplyr::group_by(country) %>%
  arrange(year) %>%
  mutate(share_ld1=lead(share,1)) %>%
  ungroup()

# Data on women in top 1 percent of the income distribution
inc <- read_excel("Womens share of top 1 income group.xlsx")
inc <- inc %>% 
  dplyr::group_by(country) %>%
  dplyr::arrange(year) %>%
  dplyr::mutate(share_t1=lead(share,1),
                share_t2=lead(share,2),
                share_t3=lead(share,3),
                share_t4=lead(share,4),
                share_avg4=(share+share_t1+share_t2+share_t3+share_t4)/5)

# Combine US and Norwegian main datasets
d <- rbind.fill(d_us,d_no)

# Add data on women in parliament
d <- merge(d,parl_full%>%
             dplyr::rename(women_share=share_avg4,
                           YEAR=year)%>%
             mutate(country=dplyr::recode(country,"us"="United States",
                                          "norway"="Norway")),
           by=c("country","YEAR"))

# Add data on share of women in top 1 percent
d <- merge(d,inc %>%
             dplyr::rename(top1_share=share_avg4,
                           YEAR=year)%>%
             mutate(country=dplyr::recode(country,"us"="United States",
                                          "norway"="Norway")),
           by=c("country","YEAR"),all.x = T)

d$outcome <- d$OUTCOME
d$all <- d$PREDALL_SW

# Calculate share of women and men who supported each policy proposal 
d$male <- d$MALE_FAV/(d$MALE_FAV+d$MALE_OPP)
d$male <- ifelse(d$SWITCHER==1,1-d$male,d$male) # Reverse the scale of the ones that should be reversed
d$male_logit <- logitTransform(d$male) # Make logit transformed version
d$female <- d$FEMALE_FAV/(d$FEMALE_FAV+d$FEMALE_OPP)
d$female <- ifelse(d$SWITCHER==1,1-d$female,d$female)
d$female_logit <- logitTransform(d$female)

d$outcome <- car::recode(d$outcome, "0=0; 1=1;2=1;3=1;4=1; else=NA")

d$inc90_logit<-logitTransform(d$PRED90_SW)
d$inc50_logit<-logitTransform(d$PRED50_SW)
d$inc10_logit<-logitTransform(d$PRED10_SW)
d$all_logit<-logitTransform(d$PREDALL_SW)

rec <- function(x) {
  x <- ifelse(d$SWITCHER==1,(1-x),x)
}

d <- d %>% mutate_at(vars(pred_i90_men:pred_i10_women),c(rec))
d <- d %>% mutate_at(vars(pred_i90_men_int:pred_i10_women_int),c(rec))

# Distance variables
d$genderdiff <- abs(d$male-d$female)
d$genderdiff_nonabs <- d$male-d$female
d$genderdiff_nonabs_i10 <- d$pred_i10_men-d$pred_i10_women
d$genderdiff_nonabs_i50 <- d$pred_i50_men-d$pred_i50_women
d$genderdiff_nonabs_i90 <- d$pred_i90_men-d$pred_i90_women

d$inc9050_diff <- abs(d$PRED90_SW-d$PRED50_SW)
d$inc9010_diff <- abs(d$PRED90_SW-d$PRED10_SW)
d$inc5010_diff <- abs(d$PRED50_SW-d$PRED10_SW)
d$inc9050_diff_nonabs <- d$PRED90_SW-d$PRED50_SW
d$inc9010_diff_nonabs <- d$PRED90_SW-d$PRED10_SW
d$inc5010_diff_nonabs <- d$PRED50_SW-d$PRED10_SW

# Mean centering variables 
d$all_mean <- d$all-mean(d$all,na.rm=T)
d$YEAR_mean <- d$YEAR-mean(d$YEAR,na.rm=T)
d$women_share_mean <- d$women_share-mean(d$women_share,na.rm=T)
d$top1_share_mean <- d$top1_share-mean(d$top1_share,na.rm=T)

# Policy area 
d$topic_comb <- ifelse(d$ECYN==1, "ep",
                       ifelse(d$SWYN==1,"ep",
                              ifelse(d$RLYN==1,"vi",
                                     ifelse(d$FPYN==1,"fp","ot"))))
d$topic_comb <- coalesce(d$topic_comb,d$TOPIC)

# Gender diff, DK included
d$male_withDK <- d$MALE_FAV/(d$MALE_FAV+d$MALE_OPP+d$MALE_DK)
d$male_withDK <- ifelse(d$SWITCHER==1,1-d$male_withDK,d$male_withDK)
d$female_withDK <- d$FEMALE_FAV/(d$FEMALE_FAV+d$FEMALE_OPP+d$FEMALE_DK)
d$female_withDK <- ifelse(d$SWITCHER==1,1-d$female_withDK,d$female_withDK)
d$genderdiff_nonabs_withDK <- d$male_withDK-d$female_withDK

d$all_withDK <- (d$MALE_FAV+d$FEMALE_FAV)/(d$FEMALE_FAV+d$FEMALE_OPP+d$FEMALE_DK+d$MALE_FAV+d$MALE_OPP+d$MALE_DK)
d$all_withDK <- ifelse(d$SWITCHER==1,1-d$all_withDK,d$all_withDK)

# DK share 
d$share_DK <- (d$MALE_DK+d$FEMALE_DK)/(d$FEMALE_FAV+d$FEMALE_OPP+d$FEMALE_DK+d$MALE_FAV+d$MALE_OPP+d$MALE_DK)
d$share_DK_male <- d$MALE_DK/(d$MALE_DK+d$MALE_FAV+d$MALE_OPP)
d$share_DK_female <- d$FEMALE_DK/(d$FEMALE_DK+d$FEMALE_FAV+d$FEMALE_OPP)


#####################################################
## Tables and Figures (same order as in appendix) ###
#####################################################


# Table A1: Gender gap in responsiveness, by policy area.
d1 <- d %>% filter(country=="United States")
d1$out1 <- d1$outcome
d1$out2 <- d1$outcome
d1$out3 <- d1$outcome
d1$out4 <- d1$outcome
m1 <- lm(out1~all+genderdiff_nonabs,data=d1%>%filter(topic_comb=="ep"))
m2 <- lm(out2~all+genderdiff_nonabs,data=d1%>%filter(topic_comb=="vi"))
m3 <- lm(out3~all+genderdiff_nonabs,data=d1%>%filter(topic_comb=="fp"))
m4 <- lm(out4~all+genderdiff_nonabs,data=d1%>%filter(topic_comb=="ot"))
d1 <- d %>% filter(country=="Norway")
d1$out5 <- d1$outcome
d1$out6 <- d1$outcome
d1$out7 <- d1$outcome
d1$out8 <- d1$outcome
m5 <- lm(out5~all+genderdiff_nonabs,data=d1%>%filter(topic_comb=="ep"))
m6 <- lm(out6~all+genderdiff_nonabs,data=d1%>%filter(topic_comb=="vi"))
m7 <- lm(out7~all+genderdiff_nonabs,data=d1%>%filter(topic_comb=="fp"))
m8 <- lm(out8~all+genderdiff_nonabs,data=d1%>%filter(topic_comb=="ot"))
stargazer(m1,m2,m3,m4,m5,m6,m7,m8, header=FALSE,
          omit.stat=c("f","ser","adj.rsq","AIC"),
          dep.var.caption = "DV: Policy adopted within 4 years of survey (1, 0)",
          dep.var.labels = c(rep("U.S.",4),rep("Norway",4)),
          column.labels=rep(rep(c("Economy","Moral","Foreign","Other")),2),
          title = "Gender gap in responsiveness, by policy area.",
          covariate.labels = c("Overall support",
                               "Diff. Men-Women"),
          column.sep.width = "-5pt",
          font.size = "small",
          notes.align = "l",
          notes.label = "",
          notes = "\\parbox[t]{1\\textwidth}{Note: Linear probability models. Variables beginning with 'Diff.' measure the support difference between two groups, calculated as (support of group A) - (support of group B). The 'Economy' area refers to policies partaining to the economy and social welfare (which are treated as seperate categories in Gilens (2012)).}")


# Table A2: How the probability of policy adoption is affected by preference difference between men and women (gender gap in responsiveness), with Dont Know’s included (support is then calculated as the number of people who favored the policy divided by the number of people who either favored, opposed, or answered dont know).
d1 <- d %>% filter(country=="United States")
m1 <- lm(outcome~all_withDK,data=d1)
m2 <- lm(outcome~all_withDK+genderdiff_nonabs_withDK,data=d1)
d1 <- d %>% filter(country=="Norway")
m3 <- lm(outcome~all_withDK,data=d1)
m4 <- lm(outcome~all_withDK+genderdiff_nonabs_withDK,data=d1)
stargazer(m1,m2,m3,m4, header=FALSE,
          omit.stat=c("f","ser","adj.rsq","AIC"),
          dep.var.caption = "DV: Policy adopted within 4 years of survey (1, 0)",
          dep.var.labels = "",
          column.labels=c("U.S.","U.S.","Norway","Norway"),
          title = "How the probability of policy adoption is affected by preference difference between men and women (gender gap in responsiveness), with Dont Know's included (support is then calculated as the number of people who favored the policy divided by the number of people who either favored, opposed, or answered dont know). ",
          covariate.labels = c("Overall support",
                               "Diff. Men-Women"),
          column.sep.width = "10pt",  
          notes.align = "l",
          notes.label = "",
          notes = "\\parbox[t]{0.75\\textwidth}{Note: Linear probability models. Variables beginning with 'Diff.' measure the support difference between two groups, calculated as (support of group A) - (support of group B).}")


# Figure B1: Policy responsiveness to men and women at different cutoff-points for preference gap
## UNITED STATES 
sequence <- seq(0.,0.12,0.01)
ests <- NULL
p <- NULL
n <- NULL
for(i in sequence){
  d1 <- d %>% filter(genderdiff>i & country=="United States")
  fit = summary(lm(outcome~male,data=d1))
  n <- rbind(n,nrow(subset(d1,!is.na(outcome))))
  ests <- data.frame(rbind(ests,fit$coefficients[,1]))
  ests$gender <- "Men"
  p <- data.frame(rbind(p,fit$coefficients[,4]))
  cutmale <- cbind(ests,p,n)[,-c(1,4)]
  names(cutmale)[1] <- "coef"
  names(cutmale)[3] <- "p"
}
sequence <- seq(0.,0.12,0.01)
ests <- NULL
p <- NULL
n <- NULL
for(i in sequence){
  d1 <- d %>% filter(genderdiff>i & country=="United States")
  fit = summary(lm(outcome~female,data=d1))
  n <- rbind(n,nrow(subset(d1,!is.na(outcome))))
  ests <- data.frame(rbind(ests,fit$coefficients[,1]))
  ests$gender <- "Women"
  p <- data.frame(rbind(p,fit$coefficients[,4]))
  cutfemale <- cbind(ests,p,n)[,-c(1,4)]
  names(cutfemale)[1] <- "coef"
  names(cutfemale)[3] <- "p"
}
cutmale$cutoff <- sequence
cutfemale$cutoff <- sequence
cuts <- rbind(cutfemale,cutmale)
cuts$sig <- car::recode(cuts$p,"lo:0.01='***'; 0.01:0.05='**'; 0.05:0.1='*'; 0.1:hi='ns'")
cuts$hjust <- ifelse(cuts$gender=="Male", 0, 1)

cuts$cutoff_n <- paste0(">",cuts$cutoff,"
                        (n=",cuts$n,")")
cuts_us <- cuts
cuts_us$country <- "United States"
## NORWAY 
sequence <- seq(0.,0.12,0.01)
ests <- NULL
p <- NULL
n <- NULL
for(i in sequence){
  d1 <- d %>% filter(genderdiff>i & country=="Norway")
  fit = summary(lm(outcome~male,data=d1))
  n <- rbind(n,nrow(subset(d1,!is.na(outcome))))
  ests <- data.frame(rbind(ests,fit$coefficients[,1]))
  ests$gender <- "Men"
  p <- data.frame(rbind(p,fit$coefficients[,4]))
  cutmale <- cbind(ests,p,n)[,-c(1,4)]
  names(cutmale)[1] <- "coef"
  names(cutmale)[3] <- "p"
}
sequence <- seq(0.,0.12,0.01)
ests <- NULL
p <- NULL
n <- NULL
for(i in sequence){
  d1 <- d %>% filter(genderdiff>i & country=="Norway")
  fit = summary(lm(outcome~female,data=d1))
  n <- rbind(n,nrow(subset(d1,!is.na(outcome))))
  ests <- data.frame(rbind(ests,fit$coefficients[,1]))
  ests$gender <- "Women"
  p <- data.frame(rbind(p,fit$coefficients[,4]))
  cutfemale <- cbind(ests,p,n)[,-c(1,4)]
  names(cutfemale)[1] <- "coef"
  names(cutfemale)[3] <- "p"
}
cutmale$cutoff <- sequence
cutfemale$cutoff <- sequence
cuts <- rbind(cutfemale,cutmale)
cuts$sig <- car::recode(cuts$p,"lo:0.01='***'; 0.01:0.05='**'; 0.05:0.1='*'; 0.1:hi='ns'")
cuts$hjust <- ifelse(cuts$gender=="Men", 0, 1)

cuts$cutoff_n <- paste0(">",cuts$cutoff,"
                        (n=",cuts$n,")")
cuts_no <- cuts
cuts_no$country <- "Norway"
## Combine and plot
cuts <- rbind(cuts_us,cuts_no)
ggplot(cuts,aes(x=reorder(cutoff_n,cutoff),
                y=coef,
                group=gender,
                shape=gender,
                color=gender)) +
  geom_line() +
  geom_point() +
  labs(x="Preference gap cutoff point",y="Effect on policy (regression coefficient)", 
       group=NULL,
       linetype=NULL,
       shape=NULL,
       color=NULL) +
  scale_color_viridis_d(begin = 0,end=0.75) +
  geom_hline(yintercept=0,color="grey60") +
  theme(legend.position = c(0.08, 0.7),
        legend.background = element_rect(size=0.2, colour ="black"),
        legend.margin = margin(1,1,1,1),
        axis.text.x =element_text(size=7.2)) +
  geom_text(aes(label=sig,hjust=hjust,color=gender),show.legend = FALSE,vjust=0) +
  facet_wrap(~country,ncol=1,scale="free_x")


############## Function for some of the next tables #################
ests <- NULL
eff<-function(x){
  m<-lm(outcome~x,data=d1)
  ests <- data.frame(rbind(ests,round(summary(m)$coefficients[2,1],2)))
  p <- summary(m)$coefficients[2,4]
  se <- round(summary(m)$coefficients[2,2],2)
  p1 <- car::recode(p,"lo:0.009999='***'; 0.01:0.049999='**'; 0.05:0.099999='*'; 0.1:hi=''")
  result <- paste0(ests, p1)
} 
eff_se<-function(x){
  m<-lm(outcome~x,data=d1)
  ests <- data.frame(rbind(ests,round(summary(m)$coefficients[2,1],2)))
  se <- round(summary(m)$coefficients[2,2],2)
  result <- paste0(ests,"\n(",se,")")
} 
eff_num<-function(x){
  m<-lm(outcome~x,data=d1)
  ests <- data.frame(rbind(ests,round(summary(m)$coefficients[2,1],2)))
  ests
} 
se<-function(x){
  m<-lm(outcome~x,data=d1)
  se <- round(summary(m)$coefficients[2,2],2)
  se
} 
p<-function(x){
  m<-lm(outcome~x,data=d1)
  p <- summary(m)$coefficients[2,4]
  substring(as.character(ifelse(p<0.001,"<<.001", round(p,3))),2)
} 
sup20<-function(x){
  m<-lm(outcome~x,data=d1)
  round(predict(m,newdata=data.frame(x=0.2),type="response"),2)
} 
sup80<-function(x){
  m<-lm(outcome~x,data=d1)
  round(predict(m,newdata=data.frame(x=0.8),type="response"),2)
} 
fact<-function(x){
  m<-lm(outcome~x,data=d1)
  diff<- predict(m,newdata=data.frame(x=0.8),type="response")/predict(m,newdata=data.frame(x=0.2),type="response")
  round(diff,1)
} 
n<-function(x){
  m<-lm(outcome~x,data=d1)
  length(m$residuals)
} 
#########################################################################


# Table B1: Policy reponsiveness to men and women when their preferences diverge.
d1 <- d %>% filter(genderdiff>0.10 & country=="United States")
tab1<-data.frame(male=t(d1 %>% summarise_at(vars(male),c(eff,se,p,sup20,sup80,fact,n))),
                 female=t(d1 %>% summarise_at(vars(female),c(eff,se,p,sup20,sup80,fact,n))))
d1 <- d %>% filter(genderdiff>0.10 & country=="Norway")
tab2<-data.frame(male=t(d1 %>% summarise_at(vars(male),c(eff,se,p,sup20,sup80,fact,n))),
                 female=t(d1 %>% summarise_at(vars(female),c(eff,se,p,sup20,sup80,fact,n))))
tab <- cbind(tab1,tab2)
row.names(tab) <- c("Effect (OLS coefficient)","Standard error","P-value", "Predicted probability of policy change at 20% support","Predicted probability of policy change at 80% support","Relative change in probability","Observations")
kable(tab,booktabs=T, align="c",linesep = '\\addlinespace', col.names = c("Men","Women","Men","Women"),
      caption = "Policy reponsiveness to men and women when their preferences diverge.") %>%
  column_spec(1,width="2.3in") %>% 
  column_spec(2:5,width="0.6in") %>% 
  row_spec(3,hline_after = T) %>%
  row_spec(6,hline_after = T) %>%
  add_header_above(c(" "=1,"United States"=2,"Norway"=2))  %>%
  footnote("Bivariate OLS models on subsetted data. Observations are the policy proposals where support diverge by more than 10 percentage points between men and women. For these subsets, the overall adoption rate is 30.9% in the US data, 35.7% in the Norwegian data. *p<0.1; **p<0.05; ***p<0.01",threeparttable = T) 


# Figure C1: Variances for income and gender from two different data sources.
var_dat <- read_excel("Variance comparison data.xlsx") %>%
  dplyr::select(year,var_sex,var_inc,source,country)
var_dat <- var_dat %>% tidyr::gather(var,val,var_inc:var_sex)
var_dat$var <- car::recode(var_dat$var,"'var_inc'='Income'; 'var_sex'='Gender';")
var_dat$source[var_dat$source=="National Election Surveys"] <- "National Election Study"
ggplot(var_dat,aes(year,val,color=var,linetype=source)) +
  geom_line() +
  facet_wrap(~country) +
  scale_color_viridis_d(begin=0,end=0.75) +
  scale_y_continuous(limits=c(0,0.26),breaks = seq(0,0.25,0.05)) +
  labs(x="Year",
       y="Variance",
       color="Variable",
       linetype="Source")


# Table C1: Full regression results for Figure 4: Policy responsiveness to men and women by income percentile and country.
d1 <- d %>% filter((inc9010_diff>0.1 | genderdiff>0.1)&country=="United States") 
tab1<-data.frame(men_10=t(d1 %>% summarise_at(vars(pred_i10_men),c(eff,se,p,sup20,sup80,fact,n))),
                 men_50=t(d1 %>% summarise_at(vars(pred_i50_men),c(eff,se,p,sup20,sup80,fact,n))),
                 men_90=t(d1 %>% summarise_at(vars(pred_i90_men),c(eff,se,p,sup20,sup80,fact,n))))
tab2<-data.frame(women_10=t(d1 %>% summarise_at(vars(pred_i10_women),c(eff,se,p,sup20,sup80,fact,n))),
                 women_50=t(d1 %>% summarise_at(vars(pred_i50_women),c(eff,se,p,sup20,sup80,fact,n))),
                 women_90=t(d1 %>% summarise_at(vars(pred_i90_women),c(eff,se,p,sup20,sup80,fact,n))))
d1 <- d %>% filter((inc9010_diff>0.1 | genderdiff>0.1)&country=="Norway") 
tab3<-data.frame(men_10=t(d1 %>% summarise_at(vars(pred_i10_men),c(eff,se,p,sup20,sup80,fact,n))),
                 men_50=t(d1 %>% summarise_at(vars(pred_i50_men),c(eff,se,p,sup20,sup80,fact,n))),
                 men_90=t(d1 %>% summarise_at(vars(pred_i90_men),c(eff,se,p,sup20,sup80,fact,n))))
tab4<-data.frame(women_10=t(d1 %>% summarise_at(vars(pred_i10_women),c(eff,se,p,sup20,sup80,fact,n))),
                 women_50=t(d1 %>% summarise_at(vars(pred_i50_women),c(eff,se,p,sup20,sup80,fact,n))),
                 women_90=t(d1 %>% summarise_at(vars(pred_i90_women),c(eff,se,p,sup20,sup80,fact,n))))
tab <- cbind(tab1,tab2,tab3,tab4)
tab1 <- data.frame(t(tab))
row.names(tab1) <- c("Income P10",
                     "Income P50",
                     "Income P90",
                     "Income P10 ",
                     "Income P50 ",
                     "Income P90 ",
                     "Income P10  ",
                     "Income P50  ",
                     "Income P90  ",
                     "Income P10   ",
                     "Income P50   ",
                     "Income P90   ")
kable(tab1[-3],booktabs=T, align="c",linesep = '\\addlinespace', col.names = c("Effect (OLS coefficient)","Standard error", "Predicted probability of policy change at 20% support","Predicted probability of policy change at 80% support","Relative change in probability","N"),
      caption = "Full regression results for Figure 4: Policy responsiveness to men and women by income percentile and country.") %>%
  column_spec(1,width="1.4in") %>%
  column_spec(2,width="0.5in") %>% 
  column_spec(3,width="0.55in") %>% 
  column_spec(4:5,width="0.8in") %>% 
  column_spec(6,width="0.5in") %>% 
  column_spec(7,width="0.4in") %>% 
  pack_rows("US: Men",1,3) %>%
  pack_rows("US: Women",4,6) %>%
  pack_rows("Norway: Men",7,9) %>%
  pack_rows("Norway: Women",10,12) %>%
  footnote("Linear probability models. Observations are the proposed policy changes where preference diverge by more than 10 percentage points either between men and women or between the 90th and 10th income percentile in each country respectively. *p<0.1; **p<0.05; ***p<0.01",threeparttable = T)


# Table C2: Policy responsiveness to men and women by income percentile and country, with support variables estimated using an interaction term
d1 <- d %>% filter((inc9010_diff>0.1 | genderdiff>0.1)&country=="United States") 
tab1<-data.frame(men_10=t(d1 %>% summarise_at(vars(pred_i10_men_int),c(eff,se,p,sup20,sup80,fact,n))),
                 men_50=t(d1 %>% summarise_at(vars(pred_i50_men_int),c(eff,se,p,sup20,sup80,fact,n))),
                 men_90=t(d1 %>% summarise_at(vars(pred_i90_men_int),c(eff,se,p,sup20,sup80,fact,n))))
tab2<-data.frame(women_10=t(d1 %>% summarise_at(vars(pred_i10_women_int),c(eff,se,p,sup20,sup80,fact,n))),
                 women_50=t(d1 %>% summarise_at(vars(pred_i50_women_int),c(eff,se,p,sup20,sup80,fact,n))),
                 women_90=t(d1 %>% summarise_at(vars(pred_i90_women_int),c(eff,se,p,sup20,sup80,fact,n))))
d1 <- d %>% filter((inc9010_diff>0.1 | genderdiff>0.1)&country=="Norway") 
tab3<-data.frame(men_10=t(d1 %>% summarise_at(vars(pred_i10_men_int),c(eff,se,p,sup20,sup80,fact,n))),
                 men_50=t(d1 %>% summarise_at(vars(pred_i50_men_int),c(eff,se,p,sup20,sup80,fact,n))),
                 men_90=t(d1 %>% summarise_at(vars(pred_i90_men_int),c(eff,se,p,sup20,sup80,fact,n))))
tab4<-data.frame(women_10=t(d1 %>% summarise_at(vars(pred_i10_women_int),c(eff,se,p,sup20,sup80,fact,n))),
                 women_50=t(d1 %>% summarise_at(vars(pred_i50_women_int),c(eff,se,p,sup20,sup80,fact,n))),
                 women_90=t(d1 %>% summarise_at(vars(pred_i90_women_int),c(eff,se,p,sup20,sup80,fact,n))))
tab <- cbind(tab1,tab2,tab3,tab4)
tab1 <- data.frame(t(tab))
row.names(tab1) <- c("Income P10",
                     "Income P50",
                     "Income P90",
                     "Income P10 ",
                     "Income P50 ",
                     "Income P90 ",
                     "Income P10  ",
                     "Income P50  ",
                     "Income P90  ",
                     "Income P10   ",
                     "Income P50   ",
                     "Income P90   ")
kable(tab1[-3],booktabs=T, align="c",linesep = '\\addlinespace', col.names = c("Effect (OLS coefficient)","Standard error", "Predicted probability of policy change at 20% support","Predicted probability of policy change at 80% support","Relative change in probability","N"),
      caption = "Policy responsiveness to men and women by income percentile and country, with support variables estimated using an interaction term") %>%
  column_spec(1,width="1.4in") %>%
  column_spec(2,width="0.5in") %>% 
  column_spec(3,width="0.55in") %>% 
  column_spec(4:5,width="0.8in") %>% 
  column_spec(6,width="0.5in") %>% 
  column_spec(7,width="0.4in") %>% 
  pack_rows("US: Men",1,3) %>%
  pack_rows("US: Women",4,6) %>%
  pack_rows("Norway: Men",7,9) %>%
  pack_rows("Norway: Women",10,12) %>%
  footnote("Linear probability models. Observations are the proposed policy changes where preference diverge by more than 10 percentage points either between men and women or between the 90th and 10th income percentile in each country respectively. *p<0.1; **p<0.05; ***p<0.01",threeparttable = T) 


# Table C3: Mean absolute support distance (in %-points) between men and women, by income group
d$abs_mw_90<-abs(d$pred_i90_men_int-d$pred_i90_women_int)
d$abs_mw_50<-abs(d$pred_i50_men_int-d$pred_i50_women_int)
d$abs_mw_10<-abs(d$pred_i10_men_int-d$pred_i10_women_int)
d1 <- d %>% filter(country=="United States")
d2 <- d %>% filter(country=="Norway")
d3<-data.frame(group=rep(c("P10 (poor)","P50 (middle)","P90 (affluent)"),2),
               country=c( rep(c("1. United States"),3), rep(c("2. Norway"),3)),
               means=round(c(mean(d1$abs_mw_10),
                             mean(d1$abs_mw_50),
                             mean(d1$abs_mw_90),
                             
                             mean(d2$abs_mw_10),
                             mean(d2$abs_mw_50),
                             mean(d2$abs_mw_90))*100,1))
d3 <- d3 %>% spread(country,means)
kable(d3,
      caption = "Mean absolute support distance (in \\%-points) between men and women, by income group",
      booktabs=T,
      linesep = '\\addlinespace',
      col.names = c("Income group","United States","Norway"),
      align = "lcc") %>%
  column_spec(1:3, width="1in") %>%
  footnote("Support for each gender/income combination imputed based on the variance/covariance procedure that included an interaction term (described in Appendix C).", threeparttable = T)


# Table C4: Policy responsiveness to men and women before and after adjusting for income
d1 <- d %>% filter( (genderdiff>0.1|inc9010_diff>0.1) & country=="United States")
tab1<-data.frame(pred_i10_men=t(d1 %>% summarise_at(vars(pred_i10_men),c(eff_num,eff,se,p,sup20,sup80,fact,n))),
                 pred_i10_women=t(d1 %>% summarise_at(vars(pred_i10_women),c(eff_num,eff,se,p,sup20,sup80,fact,n))))
d1 <- d %>% filter( (genderdiff>0.1|inc9010_diff>0.1) & country=="United States")
tab2<-data.frame(pred_i50_men=t(d1 %>% summarise_at(vars(pred_i50_men),c(eff_num,eff,se,p,sup20,sup80,fact,n))),
                 pred_i50_women=t(d1 %>% summarise_at(vars(pred_i50_women),c(eff_num,eff,se,p,sup20,sup80,fact,n))))
d1 <- d %>% filter( (genderdiff>0.1|inc9010_diff>0.1) & country=="United States")
tab3<-data.frame(pred_i90_men=t(d1 %>% summarise_at(vars(pred_i90_men),c(eff_num,eff,se,p,sup20,sup80,fact,n))),
                 pred_i90_women=t(d1 %>% summarise_at(vars(pred_i90_women),c(eff_num,eff,se,p,sup20,sup80,fact,n))))
tab_us0 <- cbind(tab1,tab2,tab3)
d1 <- d %>% filter( (genderdiff>0.1|inc9010_diff>0.1) & country=="Norway")
tab1<-data.frame(pred_i10_men=t(d1 %>% summarise_at(vars(pred_i10_men),c(eff_num,eff,se,p,sup20,sup80,fact,n))),
                 pred_i10_women=t(d1 %>% summarise_at(vars(pred_i10_women),c(eff_num,eff,se,p,sup20,sup80,fact,n))))
d1 <- d %>% filter( (genderdiff>0.1|inc9010_diff>0.1) & country=="Norway")
tab2<-data.frame(pred_i50_men=t(d1 %>% summarise_at(vars(pred_i50_men),c(eff_num,eff,se,p,sup20,sup80,fact,n))),
                 pred_i50_women=t(d1 %>% summarise_at(vars(pred_i50_women),c(eff_num,eff,se,p,sup20,sup80,fact,n))))
d1 <- d %>% filter( (genderdiff>0.1|inc9010_diff>0.1) & country=="Norway")
tab3<-data.frame(pred_i90_men=t(d1 %>% summarise_at(vars(pred_i90_men),c(eff_num,eff,se,p,sup20,sup80,fact,n))),
                 pred_i90_women=t(d1 %>% summarise_at(vars(pred_i90_women),c(eff_num,eff,se,p,sup20,sup80,fact,n))))
tab_no0 <- cbind(tab1,tab2,tab3)
d1 <- d %>% filter( (genderdiff>0.1|inc9010_diff>0.1) & country=="United States")
tab5<-data.frame(pred_men=t(d1 %>% summarise_at(vars(male),c(eff_num,eff,se,p,sup20,sup80,fact,n))),
                 pred_women=t(d1 %>% summarise_at(vars(female),c(eff_num,eff,se,p,sup20,sup80,fact,n))))
tab_us1 <- cbind(tab_us0,tab5)
d1 <- d %>% filter( (genderdiff>0.1|inc9010_diff>0.1) & country=="Norway")
tab6<-data.frame(pred_men=t(d1 %>% summarise_at(vars(male),c(eff_num,eff,se,p,sup20,sup80,fact,n))),
                 pred_women=t(d1 %>% summarise_at(vars(female),c(eff_num,eff,se,p,sup20,sup80,fact,n))))
tab_no1 <- cbind(tab_no0,tab6)
## Gap reduction in Norway when adjusting for income 
gap_no <- as.numeric(as.character(tab_no1$pred_men))[1]-as.numeric(as.character(tab_no1$pred_women))[1]
gap10_no <- as.numeric(as.character(tab_no1$pred_i10_men))[1]-as.numeric(as.character(tab_no1$pred_i10_women))[1]
gap50_no <- as.numeric(as.character(tab_no1$pred_i50_men))[1]-as.numeric(as.character(tab_no1$pred_i50_women))[1]
gap90_no <- as.numeric(as.character(tab_no1$pred_i90_men))[1]-as.numeric(as.character(tab_no1$pred_i90_women))[1]
gap_adj_no <- (gap10_no+gap50_no+gap90_no)/3
gap_red_no_pst <- (gap_no-gap_adj_no)/gap_no
gap90_no <- round(gap90_no,2)
gap_adj_no <- round(gap_adj_no,2)
## Gap reduction in United States when adjusting for income 
gap_us <- as.numeric(as.character(tab_us1$pred_men))[1]-as.numeric(as.character(tab_us1$pred_women))[1]
gap10_us <- as.numeric(as.character(tab_us1$pred_i10_men))[1]-as.numeric(as.character(tab_us1$pred_i10_women))[1]
gap50_us <- as.numeric(as.character(tab_us1$pred_i50_men))[1]-as.numeric(as.character(tab_us1$pred_i50_women))[1]
gap90_us <- as.numeric(as.character(tab_us1$pred_i90_men))[1]-as.numeric(as.character(tab_us1$pred_i90_women))[1]
gap_adj_us <- (gap10_us+gap50_us+gap90_us)/3
gap_red_us_pst <- (gap_us-gap_adj_us)/gap_us
## US
us_men <- c(as.character(tab_us1$pred_men)[2],
            as.character(tab_us1$pred_i10_men)[2],
            as.character(tab_us1$pred_i50_men)[2],
            as.character(tab_us1$pred_i90_men)[2],
            NA,NA)
us_women <- c(as.character(tab_us1$pred_women)[2],
              as.character(tab_us1$pred_i10_women)[2],
              as.character(tab_us1$pred_i50_women)[2],
              as.character(tab_us1$pred_i90_women)[2],
              NA,NA)
us_gap <- c(gap_us,gap10_us,gap50_us,gap90_us,round(gap_adj_us,2),paste0("-",round(gap_red_us_pst*100,1)," %"))
## Norway
no_men <- c(as.character(tab_no1$pred_men)[2],
            as.character(tab_no1$pred_i10_men)[2],
            as.character(tab_no1$pred_i50_men)[2],
            as.character(tab_no1$pred_i90_men)[2],
            NA,NA)
no_women <- c(as.character(tab_no1$pred_women)[2],
              as.character(tab_no1$pred_i10_women)[2],
              as.character(tab_no1$pred_i50_women)[2],
              as.character(tab_no1$pred_i90_women)[2],
              NA,NA)
no_gap <- c(gap_no,gap10_no,gap50_no,gap90_no,gap_adj_no,paste0("-",round(gap_red_no_pst*100,1)," %"))
## Together
names <- c("A. Unadjusted effects and gap",
           "B. Both at low income (P10)",
           "C. Both at median income (P50)",
           "D. Both at high income (P90)",
           "E. Average adjusted gap",
           "F. Change in gap (row A to E)")
tab <- data.frame(names,us_men,us_women,us_gap,no_men,no_women,no_gap)
kable(tab,booktabs=T,
      col.names = c("Estimate","Men","Women","Gap","Men","Women","Gap"),
      caption="Policy responsiveness to men and women before and after adjusting for income",
      linesep = '\\addlinespace') %>%
  add_header_above(c(" "=1, "United States"=3,"Norway"=3)) %>%
  row_spec(4,hline_after = T) %>%
  kable_styling(font_size = 10)


# Figure D1: Gender gap in policy responsiveness over time
d2 <- d %>% filter(country=="United States")
m1 <- lm(outcome~
           genderdiff_nonabs*YEAR+
           all_mean*YEAR,data=d2)
d2 <- d %>% filter(country=="Norway")
m4 <- lm(outcome~
           genderdiff_nonabs*YEAR+
           all_mean*YEAR,data=d2)
intmale<-interplot(m=m1,var1="genderdiff_nonabs",var2="YEAR")$data
intmale$gender <- "United States"
intfemale<-interplot(m=m4,var1="genderdiff_nonabs",var2="YEAR")$data
intfemale$gender <- "Norway"
int <- rbind(intmale,intfemale)
cols <- brewer_pal(palette="Set3")(10)
ggplot(int,aes(fake, coef1)) +
  geom_hline(yintercept=0,linetype=2,alpha=0.8,color=cols[5]) +
  geom_text(x=1968,y=-0.2,label="Equal policy responsiveness",hjust=0,color=cols[5]) +
  geom_line() +
  geom_line(aes(fake,ub)) +
  geom_line(aes(fake,lb)) +
  facet_wrap(~gender,ncol=1)+
  labs(x="Year",y="Gender gap in responsiveness
       (regression coefficient)") 


# Table D1: Full regression results for ’Table 5: Gender gap in responsiveness conditioned by time and the share of women in parliament’ for the United States.
d2 <- d %>% filter(country=="United States")
m1 <- lm(outcome~
           genderdiff_nonabs*YEAR_mean+
           all_mean*YEAR_mean,data=d2)
m2 <- lm(outcome~
           genderdiff_nonabs*women_share_mean+
           all_mean*women_share_mean,data=d2)
m3 <- lm(outcome~
           genderdiff_nonabs*YEAR_mean+
           genderdiff_nonabs*women_share_mean+           
           all_mean*YEAR_mean+
           all_mean*women_share_mean,data=d2)
m4 <- lm(outcome~
           genderdiff_nonabs*YEAR_mean+
           genderdiff_nonabs*women_share_mean+           
           all_mean*YEAR_mean+
           all_mean*women_share_mean+
           inc9050_diff_nonabs*YEAR_mean+
           inc5010_diff_nonabs*YEAR_mean+
           inc9050_diff_nonabs*women_share_mean+
           inc5010_diff_nonabs*women_share_mean,data=d2)
stargazer(m1,m2,m3,m4, header=FALSE, type='latex',
          omit.stat=c("f","ser","adj.rsq","AIC"),
          dep.var.caption = "DV: Policy adopted within 4 years of survey (1, 0)",
          dep.var.labels = "",
          covariate.labels = c("Diff. Men-Women",
                               "Year",
                               "Women MP\\%",
                               "Overall support",
                               "Diff. P90-P50",
                               "Diff. P50-P10",
                               "Diff. Men-Women x Year",
                               "Overall sup. x Year",
                               "Diff. Men-Women x W. MP\\%",  
                               "Overall sup. x W. MP\\%",
                               "Diff. P90-P50 x Year",
                               "Diff. P50-P10 x Year",
                               "Diff. P90-P50 x W. MP\\%",
                               "Diff. P50-P10 x W. MP\\%"),          
          column.labels=c("U.S.","U.S.","U.S.","U.S."),
          title = "Full regression results for 'Table 5: Gender gap in responsiveness conditioned by time and the share of women in parliament' for the United States.",
          no.space = T,
          column.sep.width = "12pt",
          font.size = "small",
          notes.align = "l",
          notes.label = "",
          notes = "\\parbox[t]{1\\textwidth}{Note: Linear probability models. Variables beginning with 'Diff' measure the support difference between two groups, calculated as (support of group A) - (support of group B). P90, P50 and P10 represent the 90th, 50th and 10th income percentile (90th is the richest). The variables overall support, year, and Women MP\\% have been centered on their overall means for the combined dataset.}")


# Table D2: Full regression results for ’Table 5: Gender gap in responsiveness conditioned by time and the share of women in parliament’ for Norway.
d2 <- d %>% filter(country=="Norway")
m1 <- lm(outcome~
           genderdiff_nonabs*YEAR_mean+
           all_mean*YEAR_mean,data=d2)
m2 <- lm(outcome~
           genderdiff_nonabs*women_share_mean+
           all_mean*women_share_mean,data=d2)
m3 <- lm(outcome~
           genderdiff_nonabs*YEAR_mean+
           genderdiff_nonabs*women_share_mean+           
           all_mean*YEAR_mean+
           all_mean*women_share_mean,data=d2)
m4 <- lm(outcome~
           genderdiff_nonabs*YEAR_mean+
           genderdiff_nonabs*women_share_mean+           
           all_mean*YEAR_mean+
           all_mean*women_share_mean+
           inc9050_diff_nonabs*YEAR_mean+
           inc5010_diff_nonabs*YEAR_mean+
           inc9050_diff_nonabs*women_share_mean+
           inc5010_diff_nonabs*women_share_mean,data=d2)
stargazer(m1,m2,m3,m4, header=FALSE, type='latex',
          omit.stat=c("f","ser","adj.rsq","AIC"),
          dep.var.caption = "DV: Policy adopted within 4 years of survey (1, 0)",
          dep.var.labels = "",
          column.labels=c("Norway","Norway","Norway","Norway"),
          covariate.labels = c("Diff. Men-Women",
                               "Year",
                               "Women MP\\%",
                               "Overall support",
                               "Diff. P90-P50",
                               "Diff. P50-P10",
                               "Diff. Men-Women x Year",
                               "Overall sup. x Year",
                               "Diff. Men-Women x W. MP\\%",  
                               "Overall sup. x W. MP\\%",
                               "Diff. P90-P50 x Year",
                               "Diff. P50-P10 x Year",
                               "Diff. P90-P50 x W. MP\\%",
                               "Diff. P50-P10 x W. MP\\%"),  
          title = "Full regression results for 'Table 5: Gender gap in responsiveness conditioned by time and the share of women in parliament' for Norway.",
          no.space = T,
          column.sep.width = "16pt",
          font.size = "small",
          notes.align = "l",
          notes.label = "",
          notes = "\\parbox[t]{1\\textwidth}{Note: Linear probability models. Variables beginning with 'Diff' measure the support difference between two groups, calculated as (support of group A) - (support of group B). P90, P50 and P10 represent the 90th, 50th and 10th income percentile (90th is the richest). The variables overall support, year, and Women MP\\% have been centered on their overall means for the combined dataset.}")


# Table D3: Gender gap in responsiveness conditioned by time and the share of women in top 1 percent of the income distribution
d2 <- d %>% filter(country=="United States")
m1 <- lm(outcome~
           genderdiff_nonabs*top1_share_mean+
           all_mean*top1_share_mean,data=d2)
stargazer(m1, header=FALSE, type='latex',
          omit.stat=c("f","ser","adj.rsq","AIC"),
          dep.var.caption = "DV: Policy adopted within 4 years of survey (1, 0)",
          dep.var.labels = "",
          column.labels=c("U.S."),
          covariate.labels = c("Diff. Men-Women",
                               "Women\\% in Top 1\\%",
                               "Overall support",
                               "Diff. Men-Women x W.\\% in Top 1\\%",
                               "Overall support x W.\\% in Top 1\\%"),          
          title = "Gender gap in responsiveness conditioned by time and the share of women in top 1 percent of the income distribution",
          no.space = T,
          font.size = "small",
          notes.align = "l",
          notes.label = "",
          notes = "\\parbox[t]{1\\textwidth}{Note: Linear probability models. Variables beginning with 'Diff' measure the support difference between two groups, calculated as (support of group A) - (support of group B). The variables overall support, and Top 1\\% have been centered on their overall means for the combined dataset. Analysis performed only for the United States since data on the share of women in the top 1 percent in Norway only go back as far as 1990.}")


# Figure D2: Gender gap in policy responsiveness as the share of women in parliament increases, non-linear effect
d2 <- d %>% filter(country=="Norway")
d2$women_share <- factor(d2$women_share)
m5 <- lm(outcome~
           genderdiff_nonabs*women_share+
           all*women_share,data=d2)
int<-interplot(m=m5,var1="genderdiff_nonabs",var2="women_share")$data
int1<-int%>%filter(fake==1)
int2<-int%>%filter(fake==0) %>% 
  slice(1)
int <- rbind(int1,int2)
int$value <- as.numeric(as.character(sub("women_share","", int$value)))
int$value[int$fake==0] <- 0.89
ggplot(int,aes(value,coef1,group=1)) +
  # geom_line(alpha=0.2) +
  geom_point(alpha=0.3) +
  geom_errorbar(aes(ymin=lb,ymax=ub),alpha=0.2) +
  geom_smooth(se=F) +
  labs(x="Share of women in parliament",y="Gender gap in responsiveness
       (regression coefficient)") +
  scale_x_continuous(limits=c(0,0.5),
                     breaks=seq(0,0.5,0.05),
                     labels=function(x)paste0(x*100,"%")) +
  coord_cartesian(ylim=c(-5,15))


# Figure E1: The United Nations Gender Inequality Index (GII), 1995-2019 for Norway, the US, and other established democracie
index <- read_excel("GII un over time.xlsx") # https://hdr.undp.org/en/indicators/68606#
index$region <- countrycode(index$Country,"country.name","ar5")
index <- index %>% filter((region=="OECD1990" | 
                             `HDI Rank`==1) &
                            `HDI Rank`!=54)
index1 <- index %>% gather(year,val,c(`1995`,
                                      `2000`,
                                      `2005`,
                                      `2010`,
                                      `2015`,
                                      `2019`))
index1$Country[index1$`HDI Rank`==1] <- "Norway"
index1$Country[index1$`HDI Rank`==17] <- "United States"
index1<-index1%>%filter(year!=2000)
ggplot() +
  geom_line(index1%>%filter(!`HDI Rank` %in% c(1,17)),
            map=aes(as.numeric(year),
                    as.numeric(val),
                    group=Country),size=0.5,color="grey85") +
  geom_line(index1%>%filter(`HDI Rank` %in% c(1,17)),
            map=aes(as.numeric(year),
                    as.numeric(val),
                    color=Country,
                    group=Country),size=0.9) +  
  scale_color_manual(values=c("#7AD151FF","#440154FF")) +
  labs(color=NULL,y="Gender Inequality Index",x="Year")
## Country ranking
#index3 <- index1 %>% group_by(Country) %>%
#  dplyr::summarise(mean=mean(as.numeric(val),na.rm=T))


# Figure E2: The V-Dem Women Political Empowerment Index, 1964-2018 for Norway, the US, and other established democracies.
index1 <- read_excel("VDEM gender power the West over time.xlsx")
index1$country_name[index1$country_name=="United States of America"] <- "United States"
index1 <- index1 %>% filter(year>1964 & year<2019)
ggplot() +
  geom_line(index1%>%filter(!country_name %in% c("Norway","United States")),
            map=aes(as.numeric(year),
                    as.numeric(v2x_gender),
                    group=country_name),color="grey85",size=0.5) +
  geom_line(index1%>%filter(country_name %in% c("Norway","United States")),
            map=aes(as.numeric(year),
                    as.numeric(v2x_gender),
                    color=country_name,
                    group=country_name),size=0.9) +
  scale_color_manual(values=c("#7AD151FF","#440154FF")) +
  labs(color=NULL,y="Women Empowerment Index",x="Year") 
## Country ranking
#index3 <- index1 %>% group_by(country_name) %>%
#  dplyr::summarise(mean=mean(v2x_gender,na.rm=T))


# Figure E3: Political interest by gender over time, in Norway
## ELECTION STUDY 1965
Valg1965 <- as.data.frame(read_dta("Valgundersøkelsen 1965.dta"))
Valg1965$income <- Valg1965$v291
Valg1965$sex <- as.numeric(car::recode(Valg1965$v264,"1=1;2=0;else=NA"))
Valg1965$int <- (3-Valg1965$v061)
Valg1965$int <- rescale(Valg1965$int,to=c(0,1))
m<-lm(int~sex,data=Valg1965)
plot1<-ggpredict(m,terms=c("sex"))
plot1$year <- "1965"
## ELECTION STUDY 2005
Valg2005 <- as.data.frame(read_dta("Valgundersøkelsen 2005.dta"))
Valg2005$income <- car::recode(Valg2005$v294,"0=NA;99999998:hi=NA")
Valg2005$sex <- as.factor(car::recode(Valg2005$v005,"1=1;2=0;else=NA"))
Valg2005$int <- Valg2005$v007
Valg2005$int[Valg2005$int==9] <- NA
Valg2005$int <- (3-Valg2005$int)
Valg2005$int <- rescale_max(Valg2005$int,to=c(0,1))
m<-lm(int~sex,data=Valg2005)
plot2<-ggpredict(m,terms="sex")
plot2$year <- "2005"
## Combine and plot
plot<-rbind(plot1,plot2)
plot$gender <- car::recode(plot$x,"1='Men';0='Women'")
ggplot(plot,aes(year,predicted,group=gender,fill=gender)) +
  geom_bar(stat="identity",position="dodge",width=0.6,color="black") +
  scale_y_continuous(limits=c(0,1)) +
  geom_errorbar(aes(ymin=conf.low,ymax=conf.high),
                position=position_dodge(width=0.6),
                width=0.1) +
  scale_fill_viridis_d(begin=0.3,end=0.95) +
  theme(legend.position = c(0.15,0.8),
        legend.background = element_rect(color="black")) +
  labs(y="Political interest (recoded to range from 0-1)",
       x="Year (Election Study Wave)",
       fill=NULL) +
  geom_text(aes(y=predicted+0.04,label=round(predicted,2)),position=position_dodge(width=1.1)) 


# Figure E4: Political interest by gender, in the United States and Norway
issp <- read_dta("ISSP Role Of Government V.dta")
issp$int <- rescale(6-as.numeric(car::recode(issp$v46,"0=NA; 8:9=NA;")),c(0,1))
issp$sex <- as.numeric(car::recode(issp$SEX, "2=0; 9=NA"))
issp$country2<-substr(sjlabelled::as_label(issp$country),4,20)
issp_no_us <- issp %>% filter(country2=="Norway" | 
                                country2=="United States")
m<-lm(int~sex*country2,data=issp_no_us)
plot<-ggpredict(m,terms=c("sex","country2"))
plot$group<-factor(plot$group,levels(factor(plot$group))[c(2,1)])
plot$sex <- car::recode(plot$x,"1='Men';0='Women'")
ggplot(plot,aes(group,predicted,group=sex,fill=sex)) +
  geom_bar(stat="identity",position="dodge",width=0.6,color="black") +
  scale_y_continuous(limits=c(0,1)) +
  geom_errorbar(aes(ymin=conf.low,ymax=conf.high),
                position=position_dodge(width=0.6),
                width=0.1) +
  scale_fill_viridis_d(begin=0.3,end=0.95) +
  theme(legend.position = c(0.15,0.8),
        legend.background = element_rect(color="black")) +
  labs(y="Political interest (recoded to range from 0-1)",
       x="Country",
       fill=NULL) +
  geom_text(aes(y=predicted+0.04,label=round(predicted,2)),position=position_dodge(width=1.1)) 




############## Function for the next tables #################
ests <- NULL
eff<-function(x){
  m<-glm(outcome~x,data=d1, family="binomial")
  ests <- data.frame(rbind(ests,round(summary(m)$coefficients[2,1],2)))
  p <- summary(m)$coefficients[2,4]
  se <- round(summary(m)$coefficients[2,2],2)
  p1 <- car::recode(p,"lo:0.009999='***'; 0.01:0.049999='**'; 0.05:0.099999='*'; 0.1:hi=''")
  result <- paste0(ests, p1)
} 
eff_se<-function(x){
  m<-glm(outcome~x,data=d1, family="binomial")
  ests <- data.frame(rbind(ests,round(summary(m)$coefficients[2,1],2)))
  se <- round(summary(m)$coefficients[2,2],2)
  result <- paste0(ests,"\n(",se,")")
} 
eff_num<-function(x){
  m<-glm(outcome~x,data=d1, family="binomial")
  ests <- data.frame(rbind(ests,round(summary(m)$coefficients[2,1],2)))
  ests
} 
se<-function(x){
  m<-glm(outcome~x,data=d1, family="binomial")
  se <- round(summary(m)$coefficients[2,2],2)
  se
} 
p<-function(x){
  m<-glm(outcome~x,data=d1, family="binomial")
  p <- summary(m)$coefficients[2,4]
  substring(as.character(ifelse(p<0.001,"<<.001", round(p,3))),2)
} 
sup20<-function(x){
  m<-glm(outcome~x,data=d1, family="binomial")
  round(predict(m,newdata=data.frame(x=0.2),type="response"),2)
} 
sup80<-function(x){
  m<-glm(outcome~x,data=d1, family="binomial")
  round(predict(m,newdata=data.frame(x=0.8),type="response"),2)
} 
fact<-function(x){
  m<-glm(outcome~x,data=d1, family="binomial")
  diff<- predict(m,newdata=data.frame(x=0.8),type="response")/predict(m,newdata=data.frame(x=0.2),type="response")
  round(diff,1)
} 
n<-function(x){
  m<-glm(outcome~x,data=d1, family="binomial")
  length(m$residuals)
} 
###############################################################



# Figure F1: Policy responsiveness to men and women when their preferences diverge by more than 10 percentage points in the US and Norway, logistic models
diff <- d %>% filter(genderdiff>0.10 & country=="United States")
m<-glm(outcome~male, data=diff, family="binomial")
newdata <- with(diff, data.frame(male=seq(0.05,0.95,by=0.01)))
p1 <- data.frame(pred=predict(m,newdata,type = "response", se.fit=T), labels= seq(0.05,0.95,by=0.01))
p1$group <- "Men"
m<-glm(outcome~female, data=diff, family="binomial")
newdata <- with(diff, data.frame(female=seq(0.05,0.95,by=0.01)))
p2 <- data.frame(pred=predict(m,newdata,type = "response", se.fit=T), labels= seq(0.05,0.95,by=0.01))
p2$group <- "Women"
p_us <- rbind(p1,p2)
p_us$group <- factor(p_us$group)
p_us$group <- factor(p_us$group,levels(p_us$group)[c(1,2)])
p_us$country <- "United States"
diff <- d %>% filter(genderdiff>0.10 & country=="Norway")
m<-glm(outcome~male, data=diff, family="binomial")
newdata <- with(diff, data.frame(male=seq(0.05,0.95,by=0.01)))
p1 <- data.frame(pred=predict(m,newdata,type = "response", se.fit=T), labels= seq(0.05,0.95,by=0.01))
p1$group <- "Men"
m<-glm(outcome~female, data=diff, family="binomial")
newdata <- with(diff, data.frame(female=seq(0.05,0.95,by=0.01)))
p2 <- data.frame(pred=predict(m,newdata,type = "response", se.fit=T), labels= seq(0.05,0.95,by=0.01))
p2$group <- "Women"
p_nor <- rbind(p1,p2)
p_nor$group <- factor(p_nor$group)
p_nor$group <- factor(p_nor$group,levels(p_nor$group)[c(1,2)])
p_nor$country <- "Norway"
plot <- rbind(p_us,p_nor)
plot$country <- factor(plot$country, levels(factor(plot$country))[c(2,1)])
ggplot(plot,aes(x=labels,y=pred.fit,color=group)) +
  geom_line(size=0.8) +
  # geom_ribbon(aes(ymin=pred.lwr,ymax=pred.upr),alpha=0.07) +
  scale_y_continuous(breaks=seq(0,1,0.2), labels=function(mean) paste0(mean*100,"%")) +
  scale_x_continuous(limits=c(0,1), breaks=seq(0,1,0.2), labels=function(labels) paste0(labels*100,"%")) +
  geom_text(data=subset(plot,labels==0.8), aes(label=group),vjust=c(-2.5,3,-2.5,3)) +
  labs(x="Percent favoring policy change",y="Predicted probability of policy change", 
       linetype=NULL, 
       color=NULL,
       title=NULL) +
  theme_bw() +
  scale_color_viridis_d(begin = 0,end=0.7) +
  coord_cartesian(ylim=c(0,0.9)) +
  theme(legend.position = "None",
        legend.background = element_rect(size=0.2, colour ="black"),
        text = element_text(size = 14),
        panel.grid = element_line(color="grey96"),
        panel.grid.minor = element_blank(),
        strip.background = element_rect(fill="grey96")) +
  guides(linetype = guide_legend(override.aes = list(size = 0.4))) +
  facet_wrap(~country)


# Table F1: Policy reponsiveness to men and women when their preferences diverge, logistic models
d1 <- d %>% filter(genderdiff>0.10 & country=="United States")
tab1<-data.frame(male=t(d1 %>% summarise_at(vars(male),c(eff,se,p,sup20,sup80,fact,n))),
                 female=t(d1 %>% summarise_at(vars(female),c(eff,se,p,sup20,sup80,fact,n))))
d1 <- d %>% filter(genderdiff>0.10 & country=="Norway")
tab2<-data.frame(male=t(d1 %>% summarise_at(vars(male),c(eff,se,p,sup20,sup80,fact,n))),
                 female=t(d1 %>% summarise_at(vars(female),c(eff,se,p,sup20,sup80,fact,n))))
tab <- cbind(tab1,tab2)
row.names(tab) <- c("Effect (logit coefficient)","Standard error","P-value", "Predicted probability of policy change at 20% support","Predicted probability of policy change at 80% support","Relative change in probability","Observations")
kable(tab,booktabs=T, align="c",linesep = '\\addlinespace', col.names = c("Men","Women","Men","Women"),
      caption = "Policy reponsiveness to men and women when their preferences diverge, logistic models.") %>%
  column_spec(1,width="2.3in") %>% 
  column_spec(2:5,width="0.6in") %>% 
  row_spec(3,hline_after = T) %>%
  row_spec(6,hline_after = T) %>%
  add_header_above(c(" "=1,"United States"=2,"Norway"=2))  %>%
  footnote("Bivariate logistic regression models on subsetted data. Included are policy proposals where support diverged by more than 10 percentage points between men and women. *p<0.1; **p<0.05; ***p<0.01",threeparttable = T) 


# Table F2: Policy reponsiveness to men and women overall, logistic models
d1 <- d %>% filter(country=="United States")
tab1<-data.frame(all=t(d1 %>% summarise_at(vars(all),c(eff,se,p,sup20,sup80,fact,n))),
                 male=t(d1 %>% summarise_at(vars(male),c(eff,se,p,sup20,sup80,fact,n))),
                 female=t(d1 %>% summarise_at(vars(female),c(eff,se,p,sup20,sup80,fact,n))))
d1 <- d %>% filter(country=="Norway")
tab2<-data.frame(all=t(d1 %>% summarise_at(vars(all),c(eff,se,p,sup20,sup80,fact,n))),
                 male=t(d1 %>% summarise_at(vars(male),c(eff,se,p,sup20,sup80,fact,n))),
                 female=t(d1 %>% summarise_at(vars(female),c(eff,se,p,sup20,sup80,fact,n))))
tab <- cbind(tab1,tab2)
row.names(tab) <- c("Effect (logit coefficient)","Standard error","P-value", "Predicted probability of policy change at 20% support","Predicted probability of policy change at 80% support","Relative change in probability","Observations")
kable(tab,booktabs=T, align="c",linesep = '\\addlinespace', col.names = c("All", "Men","Women","All", "Men","Women"),
      caption = "Policy reponsiveness to men and women overall, logistic models.") %>%
  column_spec(1,width="2.3in") %>% 
  column_spec(2:7,width="0.5in") %>% 
  row_spec(3,hline_after = T) %>%
  row_spec(6,hline_after = T) %>%
  add_header_above(c(" "=1,"United States"=3,"Norway"=3))  %>%
  footnote("Bivariate logistic regression models. The observations are proposed policy changes asked about in nationally representative surveys of the US population in the period 1964-2006 (Gilens 2012), and the Norwegian population in the period 1966-2014. The dependent variable is a dichotomous measure of whether or not the policy change was adopted within four years of the time of the survey question. *p<0.1; **p<0.05; ***p<0.01",threeparttable = T) 


# Table F3: The gender gap in responsiveness when accounting for unequal responsiveness between income groups, logistic models
d1 <- d %>% filter(country=="United States")
m1 <- glm(outcome~all+genderdiff_nonabs,data=d1, family="binomial")
m2 <- glm(outcome~all+inc9050_diff_nonabs+inc5010_diff_nonabs,data=d1, family="binomial")
m3 <- glm(outcome~all+genderdiff_nonabs+inc9050_diff_nonabs+inc5010_diff_nonabs,data=d1, family="binomial")
d1 <- d %>% filter(country=="Norway")
d1$`Norway` <- d1$outcome
m4 <- glm(outcome~all+genderdiff_nonabs,data=d1, family="binomial")
m5 <- glm(outcome~all+inc9050_diff_nonabs+inc5010_diff_nonabs,data=d1, family="binomial")
m6 <- glm(outcome~all+genderdiff_nonabs+inc9050_diff_nonabs+inc5010_diff_nonabs,data=d1, family="binomial")
stargazer(m1,m2,m3,m4,m5,m6, header=FALSE, type='latex',
          keep.stat=c("n"),
          dep.var.caption = "DV: Policy adopted within 4 years of survey (1, 0)",
          dep.var.labels = "",
          column.labels=c("US","US","US","Norway","Norway","Norway"),
          title = "The gender gap in responsiveness when accounting for unequal responsiveness between income groups, logistic models",
          covariate.labels = c("Overall support",
                               "Diff. Men-Women",
                               "Diff. P90-P50",
                               "Diff. P50-P10"),
          notes.align = "l",
          notes.label = "",
          notes = "\\parbox[t]{1\\textwidth}{Note: Logistic regression models. Variables beginning with 'Diff.' measure the support difference between two groups, calculated as (support of group A) - (support of group B). P90, P50 and P10 represent the 90th, 50th and 10th income percentile (90th is the richest).}")


#################### END  ###########################
